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Abstract — The lattice A* n is an important lattice because of its 
covering properties in low dimensions. Clarkson [1] described an 
algorithm to compute the nearest lattice point in A* n that requires 
0(n log n) arithmetic operations. In this paper, we describe a 
new algorithm. While the complexity is still 0(n log n), it is 
significantly simpler to describe and verify. In practice, we find 
that the new algorithm also runs faster. 

Index Terms — Lattice theory, nearest point algorithm, quanti- 
zation, channel coding, frequency estimation, direction-of-arrival 
estimation, synchronization 



I. Introduction 

THE study of point lattices is of great importance in 
several areas of number theory, particularly the studies of 
quadratic forms, the geometry of numbers and simultaneous 
Diophantine approximation, and also to the practical engineer- 
ing problems of quantization and channel coding. They are 
also important in studying the sphere packing problem and 
the kissing number problem [1], [2]. 

A lattice, L, is a set of points in R" such that 

L = {x e R n |x = Bw,w e Z n } 

where B is termed the generator matrix. 

The lattice A* is an interesting lattice due to its covering 
properties in low dimensions. It gives the thinnest covering in 
all dimensions up to 8 [2]. A* n has also found application in a 
number of estimation problems including period estimation 
from sparse timing data [3], frequency estimation [4] and 
direction of arrival estimation [5], 

The nearest lattice point problem is: Given y S R n and 
some lattice L whose lattice points lie in R™, find the lattice 
point xei such that the Euclidean distance between y and 
x is minimized. If the lattice is used for vector quantization 
then the nearest lattice point corresponds to the minimum 
distortion point. If the lattice is used as a code for a Gaussian 
channel, then the nearest lattice point corresponds to maximum 
likelihood decoding [6], 

Conway and Sloane [6] appear to have been the first to 
study the problem of computing the nearest lattice point in A* n . 
By decomposing A* into a union of translations of its dual 
lattice A n , they discovered an algorithm for computing the 
nearest lattice point to a given point in 0(n 2 logn) arithmetic 
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operations. Later [7], they were able to improve the execution 
time of the algorithm to 0(n 2 ) operations. 

Clarkson [1] further improved upon the work of Conway 
and Sloane and described an algorithm to compute the nearest 
lattice point that requires only 0(n\ogn) arithmetic opera- 
tions. In this paper we describe an algorithm that is similar 
to Clarkson's algorithm. Like Clarkson's algorithm, our algo- 
rithm requires O(nlogn) arithmetic operations. However, our 
algorithm and its derivation are simpler. The new algorithm, 
although of the same order of complexity, is computationally 
superior. 

We now describe how the paper is organized. Section HI] 
introduces some preliminary results and definitions. In Sec- 
tion [III] we derive all results necessary to prove that the algo- 
rithm does find the nearest lattice point. Section ITVl describes 
the algorithm. A pseudocode implementation is provided. Is 
Section [V] the arithmetic complexity of the algorithm is shown 
to be 0(n\ogn). We also tabulate some practical computation 
times that show the new algorithm to be computationally 
superior to Clarkson's original algorithm. 

II. Preliminary Theory 

Vectors and matrices are written in bold. The ith element 
in a vector is denoted by a subscript: x.- L . The transpose of a 
vector is indicated by superscript T: x T . We let 1 be a column 
vector of l's and be a column vector of zeros with a 1 in 
the ith position. 

The Voronoi region or nearest-neighbor region V(x) of a 
lattice point x is the subset of R™ such that, with respect to a 
given norm, all points in V(x) are nearer to x than to any other 
point in the lattice. The Voronoi regions are n dimensional 
polytopes [2]. 

The cubic lattice Z n is the set of n dimensional vectors with 
integer elements. The Voronoi regions of Z" are hypercubes 
of side length 1. 

The lattice A* n can be defined as the projection of the cubic 
lattice Z n+1 onto the hyperplane orthogonal to 1. This is, 



A; i = {Qx|xeZ" +1 } 
where Q is the projection matrix 

11 T 



Q 



i 



l 



(i) 



(2) 



where I is the (n+ 1) x (n + 1) identity matrix. 

Let II be a permutation matrix. Observe the following 
elementary properties: 

i) m = i, 
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2) i T n = i T , 

3) ||nx|| = |]x]|. 

Lemma 1: The matrices II and Q commute, i.e., II Q 
QII 

Proof: Using the properties of the permutation matrix, 
observe that 

n Q = n(i i^) = n nilT 



where we set 



n 

(i 



n + l' 
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n + 1 n+l 
11 T 

— — )n = qn. 
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Corollary 1: For all z e R"+\ ||Qz|| = ||QIIz||. 

Corollary 2: x 6 A* n if and only if IIx € A*. 

Proof: Because the inverse of a permutation matrix is 
also a permutation matrix, we need only prove sufficiency. If 
x e A* n then x = Qk with k e Therefore, IIx = 

IIQk = Qllk = Qk' where k' = nk e Z n+1 and so 

nx e A*. ■ 

Corollary 3: The lattice point x is a closest point in A* n to 
y if and only if IIx is a closest point in A* to Tly. 

Proof: As for Corollary [2] we need only show sufficiency. 
We do this by contradiction. Suppose Tlx is not closest to Ily 
but there is instead some IIz 6 A* such that 



|n(z-y)|| < ||n(x-y) 



This implies that 



y < ll x - 



which contradicts the assumption that x is a closest point to 
yinA* n . ' ■ 

Hence, in considering an algorithm to find a closest point 
in A* n to y, it is sufficient to consider a canonical permutation 
of y. We will see that it is very convenient to consider the 
permutation in which the {centered) fractional parts ofy, i.e., 
{Hi} = Vi ~ \_Vi\> are sorted in descending order. That is, in 
the sequel, except where otherwise noted, we will assume that 



{yi} > {2/2} >•••> {y n+ i}. 



(3) 



In the case that two or more {y{\ are equal then multiple 
orderings of y satisfy (|3). The following arguments and the 
subsequent algorithm are valid for any ordering of y that 
satisfies (01. 

III. Closest Point in A* n 

Lemma 2: If x = Qk is a closest point in A* n to y 6 R ra+1 
then there exists some A € K for which k is a closest point 
in Z" +1 to y + Al. 

Proof: Decompose y into orthogonal components Qy and 
tl for some t£l. Then 



Observe that 



|Q(y-k)|| 2 +t 2 (n + l). 



Q(y-k) =y + Al-k 



(4) 



A = 



l T (k-y) 
n + l 



Suppose k is not a closest point in Z™ +1 to y + Al. Suppose 
k' is closer. Let x' = Qk'. Then 

||y-x'|| 2 = ||Q(y-k')|| 2 + i 2 (n + l) 
< ||y + Al -k'|| 2 + t 2 (n+l) 
<||y + Al-k|| 2 + * 2 (n + l) = ||y-x|| 2 , 

contradicting the assumption that x is a closest point in A* n 
to y. ■ 

Now consider the function f : M 1— > Z" +1 defined so that 

f (A) = Ly + All 

where [•] applied to a vector denotes the vector in which each 
element is rounded to a nearest integeiQ. That is, f (A) gives a 
nearest point in to y + Al as a function of A. Observe 

that f (A + 1) = f (A) + 1. Hence, 



Qf(A + 1) = Qf(A). 



(5) 



Lemma [2] implies there exists some A £ K such that x = 
Qf(A) is a closest point to y. Furthermore, we see from © 
that A can be found within an interval of length 1 . Hence, if 
we define the set 

^ = {f(A)|Ae[0,l)} 

then Q^ contains a closest point in A* n to y. 

If the fractional parts of y are sorted as in y), it is clear 
that 5? contains at most n + 2 vectors, i.e., 

^C{[yl, Lyl+ei, Lyl+ei+e 2 ,...,Lyl+e 1 +---+e n+ i}. 

(6) 

It can be seen that the last vector listed in the set is simply 
Ly] + 1 and so, once multiplied by Q, the first and the last 
vector are identical. 

An algorithm immediately suggests itself: test each of the 
n + 1 distinct vectors and find the closest one to y. Indeed, 
this is exactly the principle of the algorithm we propose here. 
It only remains to show that this can be done in 0(n log n) 
arithmetic operations. 

IV. Algorithm 

We label the elements of 5^ according to the order given 
in ([SJ. That is, we set u = |_y] and, for i = 1, . . . , n, 



u. t = u. t _i + e I 



(7) 



Let Zj = y — U;. Clearly, zq = {y}. Following (0]), the squared 
distance between Qu; and y is 



Hy-QuiH 2 =di+t 2 (n + l) 
where we define d, as 



di = ||Qzj 



z/1 ■ 

n + l' 



Mi)! 

n+l 



= z, %i - 



(8) 



(9) 



' The direction of rounding for half-integers is not important. However, the 
authors have chosen to round up half-integers in their own implementation. 
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We know that the nearest point to y is that Qiij which 
minimizes ((H). Since the term t 2 (n + 1) is independent of the 
index i, we can ignore it. That is, it is sufficient to minimize 
di, i = 0, . . . , n. 

We now show that di can be calculated inexpensively in a 
recursive fashion. We define two new quantities, on = zf 1 
and fii = zjzi. From ((7), 

a. t = zf 1 = ( Zi _! - e;) T l = ai _ a - 1 (10) 

and 

A = zf Zj = (zj_ x - e. i ) T (z i _i - e 4 ) = - 2{y.;} + 1. 

(11) 



This appears to be the case for small n. Table U shows the 
practical computational performance of Clarkson's algorithm 
versus our new algorithm. It is evident that the new algorithm 
is computationally superior, particularly for small n. It appears 
that the computational performance of the algorithms converge 
for large n. The computer used for these trials is an Intel Core2 
running at 2.13Ghz. 

TABLE I 

Computation time in seconds for 10 s trials 

Algorithm n=20 n=50 n=100 n=500 
Clarkson 457 657 11.11 47.81 
New 2.05 3.86 7.125 35.44 



Input: y G R n+1 

1 z = y - iyi 

2 a = z T l 

3 (3 = z T z 

4 s = dsortindices(z) 

6 m = 

7 for i = 1 to n do 

a = a — 1 
P = P- 2z Si 4 



9 
10 

11 

12 



1 



if (3 



n+l 

D = f3 
m = i 



< D then 



n+l 



13 k = Lyl 

14 for i = 1 to m do 

is I fc Ss = k Si + 1 



16 x = k 

17 return x 



i J k- 



Algorithm 1: Algorithm to find a nearest lattice point in A„ to 

y G R n+1 



Algorithm [T] now follows. The main loop beginning at 
line |7] calculates the a, and /3j recursively. There is no 
need to retain their previous values, so the subscripts are 
dropped. The variable D maintains the minimum value of the 
(implicitly calculated values of) di so far encountered, and m 
the corresponding index. 



As a final note, the algorithm proposed here can be extended 
to other lattices for which Lemmata [T] and |2] hold. Potential 
candidates are the Coxeter lattices [8], [9]. 
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V. Computational Complexity 

Each line of the main loop requires O(l) arithmetic compu- 
tations so the loop (and that on line (14) requires O(n) in total. 
On line|4]the function dsortindices(z) returns the vector s such 
that z Sl > z S2 > ■ ■ ■ > z Sn+1 . This sorting operation requires 
0(n log n) arithmetic operations. The vector operations on 
lines [lfO Q3] and [16] all require 0(n) operations. It can be 
seen, then, that the computational cost of the algorithm is dom- 
inated by the sorting operation and is therefore O(nlogn). 

Clarkson's original algorithm required two sorts of n + 1 
elements. The new algorithm requires only a single sort. See- 
ing as the sort dominates the complexity of both algorithms, 
we might expect our algorithm to require approximately half 
the arithmetic operations of Clarkson's original algorithm. 



